#!/usr/bin/perl
use Cwd 'abs_path';
use Getopt::Std;
#Input Arguments
getopt("d:t:a:f:o", \%args);
#my $dir_or = $args{d};

my $tped= $args{t};
my $arraytype= $args{a};
my $tfam = $args{f};
my $outdir = $args{o};
my $dirtemp = $args{d};
chomp($dirtemp);
if($dirtemp eq "")
{
        die "entered temp dir name  is empty\n";

}
if (!(-e $tfam)) {
        die "input tfam file does not exist\n";
}
if (!(-e $tped)) {
        die "input tped file does not exist\n";
}
if($tped eq "")
{
        die "entered tped file name  is empty\n";

}
else
{
@tped_ar=split(/\//,$tped);
$tped = pop(@tped_ar);
$dir_or = join('/',@tped_ar);
}
if($tfam eq "")
{
	die "entered tfam file name  is empty\n";
}
if($outdir eq "")
{
	die "entered output directory name  is empty\n";
}

#to know the path of this script to grab working directory location
$line = abs_path($0);
chomp $line;

@DR_array = split('/',$line);
pop(@DR_array);
$dir = join("/",@DR_array);
$round = sprintf("%.0f", rand()*time());
$rounded = "temp".$round;
#$rounded = "temp317759169";
#creating a temp directory if not presen
system("mkdir $dirtemp");
system("mkdir $dirtemp/$rounded");
system("cp $tfam  $dirtemp/$rounded");
@tfam_ar=split(/\//,$tfam);
$tfam = pop(@tfam_ar);
system("mv $dirtemp/$rounded/$tfam $dirtemp/$rounded/unprocessed_beagle_input.tfam");
system("cp $dir_or/$tped $dirtemp/$rounded/unprocessed_beagle_input.tped");
if($arraytype =~ m/AFFY/)
{
	$geno = 0.1;
	$mind = 0.1;
}
else
{
	$geno = 0.05;
	$mind = 0.05;
}
print "geno: $geno mind: $mind\n";

$tfam ="unprocessed_beagle_input.tfam";
$tped = "unprocessed_beagle_input.tped";
system("$dir/bin/plink/plink  --tfile $dirtemp/$rounded/unprocessed_beagle_input --missing --out $dirtemp/$rounded/miss");
if (!(-e "$dirtemp/$rounded/miss.imiss")) {
        die "something wrong in the plink command miss \n";
}
system("$dir/bin/plink/plink  --tfile $dirtemp/$rounded/unprocessed_beagle_input --missing --geno 0.999 --mind 0.15 --make-bed --out $dirtemp/$rounded/intialqc");
if (!(-e "$dirtemp/$rounded/intialqc.fam")) {
        die "something wrong in the plink command geno\n";
}

system("$dir/bin/plink/plink  --bfile $dirtemp/$rounded/intialqc --missing --geno $geno --make-bed --out $dirtemp/$rounded/geno");
if (!(-e "$dirtemp/$rounded/intialqc.fam")) {
        die "something wrong in the plink command geno\n";
}

system("$dir/bin/plink/plink --bfile $dirtemp/$rounded/geno --missing --mind $mind --make-bed --out $dirtemp/$rounded/mind");
if (!(-e "$dirtemp/$rounded/mind.fam")) {
        die "something wrong in the plink command mind\n";
}

system("$dir/bin/plink/plink  --bfile $dirtemp/$rounded/mind --impute-sex  --recode --transpose --out $dirtemp/$rounded/imputesex");
if (!(-e "$dirtemp/$rounded/imputesex.tfam")) {
        die "something wrong in the plink command impute sex\n";
}
system("mv $dirtemp/$rounded/miss.imiss $dirtemp/$rounded/samplesmiss.qc");
#system("mv $dir/temp/$rounded/samples.qc $outdir/samples.qc");
system("rm $dirtemp/$rounded/miss.*");
system("rm $dirtemp/$rounded/geno.*");
system("rm $dirtemp/$rounded/mind.*");
system("mv $dirtemp/$rounded/imputesex.sexcheck $dirtemp/$rounded/sexcheck.txt");
system("mv $dirtemp/$rounded/imputesex.tped $dirtemp/$rounded/unprocessed_beagle_input.tped");
$sys = "perl $dir/bin/perl_replace_originalsex_predicted_sex.pl $dirtemp/$rounded/sexcheck.txt $dirtemp/$rounded/imputesex.tfam $dirtemp/$rounded/unprocessed_beagle_input.tfam";
print $sys."\n";
system($sys);
system("rm $dirtemp/$rounded/imputesex.*");
system("cp $dirtemp/$rounded/unprocessed_beagle_input.tfam $dirtemp/$rounded/processed_beagle_input.tfam");
$sys = "perl $dir/bin/sexcheck.pl -i $dirtemp/$rounded/samplesmiss.qc -d $dirtemp/$rounded/sexcheck.txt -v $mind -o $dirtemp/$rounded/samples.qc";
print $sys."\n";
system($sys);
#system("mv $dir/temp/$rounded/sexcheck.txt $outdir/sexcheck.txt");
system("mv $dirtemp/$rounded/samples.qc $outdir/sample.qc");
=head
#system("mv $dir/temp/$rounded/genotypes.qc $outdir/genotypes.qc");
#forward strand mapping
@array_arraytype = split(',',$arraytype);
if(@array_arraytype > 1)
{
	$refpath = "$dir/bin/FwdStrand_Naresh/ref/".$array_arraytype[0]."/FWD_REF_INFO";
	for($ij=1;$ij<@array_arraytype;$ij++)
	{
		$refpath =$refpath.","."$dir/bin/FwdStrand_Naresh/ref/".$array_arraytype[$ij]."/FWD_REF_INFO";
	}
}
else
{
	$refpath = "$dir/bin/FwdStrand_Naresh/ref/$arraytype/FWD_REF_INFO";
}
#$cmd_system = "/bin/sh $dir/bin/FwdStrand_Naresh/bin/RunFwd.sh $dir/bin/FwdStrand_Naresh $dir/temp/$rounded/ $dir_or/$tped $dir/temp/$rounded/logfile $arraytype $refpath";
$cmd_system = "/bin/sh $dir/bin/FwdStrand_Naresh/bin/RunFwd.sh $dir/bin/FwdStrand_Naresh $dirtemp/$rounded/ $dirtemp/$rounded/unprocessed_beagle_input.tped $dirtemp/$rounded/logfile $arraytype $refpath";
print $cmd_system."\n";
system($cmd_system);
wait;
#system("cp $tfam  $dir/temp/$rounded");
#@tfam_ar=split(/\//,$tfam);
#$tfam = pop(@tfam_ar);

wait;
$fwdstrd_output_tped = "fwdStrandResults_$tped";
$fwdstrd_output_ind = "fwdStrandFlippedInd_$tped";
$fwdstrd_output_ind =~ s/.tped/.txt/g ;
#separating the snps which are ambigous and snps that are not present in the hapmap3 based on fwd strand results
wait;

$command_sys = "/usr/bin/perl $dir/bin/process_fwdstnd_beagle_input.pl -f $dirtemp/$rounded/$fwdstrd_output_tped -t $dirtemp/$rounded/$fwdstrd_output_ind -h $dir/bin/files/newhapmap3_sample_run.tped  -n $dirtemp/$rounded/processed_beagle_input.tped -e $dirtemp/$rounded/excluded_no_hapmap_processed_beagle_input.tped -a $dirtemp/$rounded/ambi_hapmap_processed_beagle_input.tped -u $dirtemp/$rounded/unsure_fwdmapping.tped";
print $command_sys."\n";
system($command_sys);
#system("mv $dir/temp/$rounded/$tfam $dir/temp/$rounded/processed_beagle_input.tfam");
#system("cp $dir_or/$tped $dir/temp/$rounded/unprocessed_beagle_input.tped");
#creating binary plink files
system("$dir/bin/plink/plink --tfile $dirtemp/$rounded/processed_beagle_input --out $dirtemp/$rounded/processed_beagle_input --make-bed");
#chopping the chr 1-22 and submitting the jobs for imputation
$syst = "/usr/bin/perl $dir/bin/perl_create_ind_chr_imputation_ori.pl -d $dir -t $dirtemp/$rounded";

print $syst."\n";
system($syst); 
###$sys = "/usr/bin/perl $dir/bin/perl_create_ind_chr_imputation_ori_X.pl -d $dir -t $rounded";
#submitting the chr x
####system($sys);
open(BU,"$dirtemp/$rounded/submitjobs");
open(WRBU,"$dirtemp/$rounded/submitjobs1");
$line_jobs = "";
$submitjobs_num = 0;
while($line = <BU>)
{
	chomp($line);
	$submitjobs_num++;
	$line = $line." >>$dirtemp/$rounded/submitjobs1";
	system($line);
	
}
$submitjobs_num1 = 0;
open(BU1,"$dirtemp/$rounded/submitjobs1");
while($line = <BU1>)
{
        chomp($line);
	$submitjobs_num1++;
	@array = split(" ",$line);
	$line_jobs = $line_jobs.$array[2]." ";
}
if($submitjobs_num != $submitjobs_num1)
{
	die "all jobs are not submitted\n";
}

#checking the submitted jobs which are completed
system("/home/inforsense/applications/groovy-1.6.4/bin/groovy  $dir/bin/CheckJobStatus.groovy $line_jobs");
#merging the result files chromosome wise 
#open(WWW, ">$dic/temp/$rounded/failed.txt");
open(WWW, ">$dirtemp/$rounded/failed.txt");
$sys= "/usr/bin/perl $dir/bin/perl_after_imputation_merge.pl -d $dir -t $dirtemp/$rounded";
system($sys);
#creating the combined male and female  x chr imputation results files
###system("/usr/bin/perl $dir/bin/perl_club_male_female_x_chr_after_imp.pl -d $dir/temp/$rounded");
open BUF,"< $dirtemp/$rounded/failed.txt" or die "can't open failed.txt";
#checking if any chr outfile is not created or failed
while($lined = <BUF>)
{
	if($lined =~ m/failed/)
	{
		system("rm $dirtemp/$rounded/whole_geno_phased");
		system("rm $dirtemp/$rounded/whole_gprobs");
		system("rm $dirtemp/$rounded/whole_no_snp_hapmap.txt");
		system("rm $dirtemp/$rounded/whole_r2");
		die "$lined\n";
	}
}
#$sys= "/usr/bin/perl $dir/bin/perl_after_imputation_merge.pl -d $dir";
#system($sys);
#system("/usr/bin/perl $dir/bin/perl_club_male_female_x_chr_after_imp.pl -d $dir/temp/$rounded");
#combing imputation output, checking ambigous snps,creating dosage and indicator file
###system("/usr/bin/perl $dir/bin/perl_club_male_female_after_imp_b4_imp.pl -d $dir  -t $rounded");
#restrict the output to the hapmap3 or 3 datasets
#system("/usr/bin/perl $dir/bin/perl_restrict_allfiles_hapmap3_after_impute.pl -t $dir/temp/whole_geno_1_22_x_nosnphapmap_nohapmapsnpsInital_check_ambi1.tped -d $dir/temp/whole_chr_gprobs_dosage_ambi -i $dir/temp/whole_ind -r $dir/files/rsids_pancreas_emerge_vte  -h $dir/files/snps_hapmap_x -o $outdir");
#system("cp $dir/temp/processed_beagle_input.tfam $outdir/final.tfam");
$maindir = "$dirtemp/$rounded";
system("/usr/bin/perl $dir/bin/perl_create_beagle_output_plink_tped.pl -m $dir/bin/files/newhapmap3_sample_run.tped -i $maindir/whole_geno_phased -o $maindir/whole_geno_phased.tped");
system("/usr/bin/perl $dir/bin/perl_create_combine_dosage_1stphase_nox.pl $maindir/whole_gprobs $maindir/whole_gprobs1 $maindir/whole_geno_phased.tped");
system("/usr/bin/perl $dir/bin/perl_create_indicator_file_1stphase_nox.pl $maindir $maindir/whole_r2 $maindir/unprocessed_beagle_input.tped $maindir/whole_geno_phased.tped $maindir/whole_gprobs1");

$syste = "/usr/bin/perl $dir/bin/perl_create_new_map_tpedfile_4sort.pl $maindir/whole_geno_phased.tped $maindir/snp_order.txt $maindir";
#print $syste."\n";
system($syste);
system("/usr/bin/perl $dir/bin/perl_inputlist_sortfile.pl $maindir/snp_order.txt $maindir/whole_gprobs1 1 $maindir/sorted_whole_gprobs1");
system("/usr/bin/perl $dir/bin/perl_inputlist_sortfile.pl $maindir/snp_order.txt $maindir/whole_ind  1 $maindir/sorted_whole_ind");
system("/usr/bin/perl $dir/bin/perl_inputlist_sortfile_tped.pl $maindir/snp_order.txt  $maindir/whole_geno_phased.tped 2 $maindir/sorted_whole_geno.map");
system("mkdir -p $outdir");
system("mv $maindir/sorted_whole_gprobs1 $outdir/final.dosage");
system("mv $maindir/sorted_whole_ind $outdir/final.ind");
system("mv $maindir/sorted_whole_geno.map  $outdir/final.map");
system("mv $maindir/processed_beagle_input.tfam $outdir/final.tfam");
#system("rm -rf $maindir/*");
#system("rm -rf $maindir");
